Determining the Spectral Signature of Spatial Coherent Structures 
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We applied to an open flow a proper orthogonal decomposition (pod) technique, on 2D snapshots 
of the instantaneous velocity field, to reveal the spatial coherent structures responsible of the self- 
sustained oscillations observed in the spectral distribution of time series. We applied the technique 
to 2D planes out of 3D direct numerical simulations on an open cavity flow. The process can easily 
be implemented on usual personal computers, and might bring deep insights on the relation between 
spatial events and temporal signature in (both numerical or experimental) open flows. 
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One of the most challenging questions arising in open 
flows such as jets, mixing layers, etc, is to understand the 
occurrence and nature of robust and reproducible self- 
sustained oscillations revealed in spatially localized time 
series, usually velocity or pressure measurements. How 
such frequencies appear, and whether or not they might 
be the signature of particular coherent spatial patterns, 
still remain largely unresolved, although abundantly doc- 
umented [ll 1^1 • Such an understanding may moreover ap- 
pear of the upmost importance in control applications, in 
that knowing which spatial event is generating such spec- 
tral signature may lead to best fitted control scheme with 
respect to the required goal. An example is given by flows 
over open cavities, like in high speed trains, that generate 
very powerful self-sustained oscillations that appear to be 
the main source of noise emitted by the train. In that 
case, control will be aimed to reduce or even suppress 
the source of noise, without reducing the aerodynamic 
performances, and at the lowest energetic cost. 

In this paper we (i) show in a test case the ability 
of the pod technique to associate self-sustained oscilla- 
tions to well- identified spatial coherent structures; (ii) 
confirm, as a consequence, the mixing layer origin of the 
most energetic self-sustained oscillations in an open cav- 
ity flow. We will show that 2D cuts out of the fully 
3D flow are sufficient to extract significant space-time 
events out of the flow. We are using for that purpose 
a technique based on an empirical decomposition of the 
flow, that optimizes a basis of (orthogonal) eigen-modcs 
with respect to the kinetic energy. The technique is of- 
ten known as to the proper orthogonal decomposition 
{pod hereafter) in the framework of fluid dynamics 
or as the Karhunen-Loeve decomposition in the frame- 
work of signal processing Q. (Other denominations ex- 
ist, such as empirical orthogonal decomposition, singular 
value decomposition, etc, depending on the field of ap- 
plication considered). To illustrate our point, we applied 
the technique to 3D direct numerical simulations of an 
air flow over an open cavity The system is a cavity 
of length L = 10 cm along x (the longitudinal direction 
along which air is flowing) , of depth h = 5 cm (the aspect 
ratio L/h is 2), and transverse size I — 20 cm. The cavity 



is enclosed into a vein 12 cm high. The flow rate velocity 
is Uq — 1.2 m/s (Reynold's number Re ~ 8500). Simula- 
tions were performed following a finite volume approach 
under an incompressible flow hypothesis. Spatial and 
time discretization have a second order precision. The 
pressure field is given by a Poisson's equation that re- 
quires a projection step, such as to be in agreement with 
a non divergent velocity field. In order to reduce the CPU 
time cost, the spanwise boundary conditions are periodic. 
The 256 x 128 x 128 mesh-spatial grid is refined in areas 
featuring strong velocity gradients (boundary and shear 
layers) — with a mesh varying from 0.7 and 10 mm along 
the longitudinal x and vertical y directions, and constant 
with about 1.56 mm over the transverse direction z [6(. 

Here we briefly expose the pod technique we im- 
plemented. The goal is to compute the eigenmodes 
{<fi n (t) , ip n (r)} that best fit the coherent structures com- 
posing the flow, computed from a data base of M dif- 
ferent snapshots of the velocity field, in such a way that 
any instantaneous snapshot of the data base can be re- 
constructed by performing the sum over the eigenmode 
basis: 
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where the A„ = //^ are the eigenvalues of the decomposi- 
tion Typically M was of the order of 600 frames. 
Note that u being a vector field, ip must also be so; 
however we will also use the notation ip when dealing 
with one component of the field (usually it will be the 
longitudinal component along x). A coherent structure 
can now be defined as an eigenmode of a (2-pointwisc 
linear) correlation matrix built on the data base snap- 
shots. There exists mainly two ways of building up a 
correlation matrix: either performing a time correlation, 
or a space correlation. With snapshots u(r, t) of size 
N = N x x N y pixels (where N x ~ 125 and N y ~ 100 are 
respectively the snapshot dimensions along x and y) , the 
space-correlation matrix 

K(r,r') = / u p (f,t)u q (f',t)dt 



is of size 2N 2 (u Ptq are velocity components). We re- 
stricted our analysis to the x,y-components of the veloc- 
ity field so as to mimic what is available from 2D ex- 
perimental PIV snapshots. On the contrary, the time 
correlation matrix 

C(t,t') = J J u{r,t)u(?,t')dr 

is of size M 2 (if M is the number of instant under con- 
sideration), much smaller than (27V) 2 (3.6 x 10 4 against 
4 x 10 s ). Keeping in mind that no more information can 
be extracted from that contained in the data base itself, 
and that at most M relevant eigenmodes are therefore 
available from the data set, we chose the second way 
(based on C(t, i')), known as the snapshot pod technique 
in the literature 0, • Practically, we start with a data 
base of M instantaneous spatial snapshots of the velocity 
field; in experiments they can for example be obtained us- 
ing PIV techniques • Then, the data are reshaped into 
a "data matrix" A whose column elements are the pixels 
of a given snapshot. For that purpose, each 2D snapshot 
is reshaped into a column vector (of length N) , by stack- 
ing over each other all the columns of the snapshot, from 
the first to the last. Both x and y components of the (vec- 
tor) velocity field are further stacked in the same column 
following the same procedure, starting with component 
x at the top of the column, and then the component y 
down to the bottom of the column. The vertical size 
of A is therefore 2N. The matrix A contains as many 
columns as snapshots in the data base (so that its hori- 
zontal dimension is M), the snapshots being ranked from 
the left to the right of A as the time is flowing down. The 
matrix A is therefore of dimension M x 2N. The correla- 
tion matrix C is next obtained by performing the product 
C — A 1 ■ A, where A* is the transposed matrix of A, and 
• the usual matrix dot product. (Note that the space 
correlation matrix K is given by K = A ■ A 1 ) Applying 
a singular value decomposition procedure on C, we ob- 
tain the eigen- modes 4> n {t), rearranged as columns of a 
chronos |9j matrix $ from left with n = 1 to right with 
n = M. The spatial eigenmodes 4>n{f) (sometimes called 
topos in the literature Q) are given following Eq.QJ by 
i'ni'r) = jj- f (f> n (t)u(r,t)dt. The ip n are reshaped into 
columns of a topos matrix ^ — (A ■ <£>) • D -1 / 2 , follow- 
ing the same procedure as A, where D is the diagonal 
matrix of the eigenvalues A n , ranked from the largest to 
the smallest value. The Matlab© software is dedicated 
to matrix operations, so that the whole process of build- 
ing A, calculating C, performing the singular decompo- 
sition to obtain the </>„, and determining the ip n , takes, 
for M = 600 and N ~ 37300 no more than 30 sec on a 
usual PC. 

We first present in FigQ the spectral distribution of 
time series provided by local recordings of one compo- 
nent of the velocity field (here the longitudinal compo- 
nent u x (t)). Velocity recordings are done at 4 different 
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FIG. 1: Power spectral distribution of the ^-component ve- 
locity time series collected in the mixing layer, upstream and 
downstream (top left to right), within the cavity, upstream 
and downstream (bottom left to right), from 3D direct nu- 
merical simulations. 



locations: 2 within the mixing layer (one upstream, one 
downstream), and 2 within the cavity (upstream and 
downstream). In each of them clearly appear peaks at 
about /o = 13.5 Hz (Strouhal number St=1.06 when 
based on the cavity length L and the reference velocity; 
St=0.033 when based on the mixing layer thickness and 
the mean velocity — to be compared with the natural 
Strouhal number St„=0.03 of an unforced mixing layer 
0), and it is now well accepted that this frequency is 
produced by the instability of the mixing layer . The 
spectral component is recovered anywhere in the cavity, 
presumably due to the overall pressure field coupling due 
to the fluid incompressibility (the Mach number is about 
4 x 10~ 3 ). 

Now we propose to apply our technique so as to iden- 
tify the spatial coherent structures ifi n (r) of the flow, and 
track out their dynamical features from their associated 
time-dependent amplitudes <j) n {t). Note that the snap- 
shots here must be sampled at least at 2/o ~ 30 Hz if 
we want to be time-resolved with respect to fo (Shannon 
criterion). This was actually achieved in the numerical 
simulations. 

In FigElwe clearly see that the pod decomposes the 
flow into two well-defined areas: one is the mixing layer 
over the cavity, essentially captured by the 2 first eigen- 
modes ipi,2, the other is the cavity vortices, captured 
by the higher order (less energetic) eigenmodes. The 
two first modes look very similar, and actually could be 
phase squared as expected when the flow is experiencing 
a global mean advection (phase squaring resulting in that 
case from the space translation invariance) 0. However, 
when comparing the eigenvalues Ai,2 plotted in FigEjfc, 
they appear to be rather different, and not close to each 
other as it should be expected in a phase squaring sit- 
uation. Moreover, when plotting chronos 4>2{t) vs 4>\{t) 
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FIG. 2: Six first spatial eigenmodes (topos tpn(r) with n — 1 
to 6 from top left to bottom right). Arrows represent the 
velocity vector in the plane of the mode (here components x 
and y). 
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FIG. 3: Singular value decomposition spectrum. Note the two 
first eigenvalues, that are very different (left). Phase portrait 
of 0a (t) versus 4>i(t) (right). 



(FigEb) , a torus is drawn whose dispersion cannot be ex- 
plained by numerical noise. Henceforth, it rather looks 
like if the two first modes were not two degenerated phase 
aspects of a unique "complex" mode, but really two dif- 
ferent pod modes, although somehow coupled (so as to 
produce the torus shape of FigEb). This invoques a sym- 
metry breaking in the flow advection, most likely due to 
the downstream corner of the cavity. In the discussion on 
whether the instability is convective or absolute, note the 
downstream corner location of the two first topos tpi,2, 
whose amplitude is vanishing in the upstream area. This 
is a strong argument in favor of the convective nature 
of the instability, the upstream front of the instability 
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FIG. 4: Five first time eigenmodes (chronos <j> n (t) with n = 1 
to 5 from top to bottom), and their associated power spectral 
distribution m. 



wavepacket being expected to spread back against the 
flow advection in an absolutely unstable situation. A 
global mode cannot be completely excluded however 0]. 

In Fig0]are shown the five first time series <j> n (t) and 
their spectral distribution. We clearly see the occurrence 
of the frequency /o = 13.5 Hz associated with the two 
first chronos, the corresponding topos featuring the co- 
herent structures contained in the mixing layer. Clearly, 
the frequency turns out to be associated with the insta- 
bility that develops in the mixing layer. While local time 
series all produce spectral components at fo (see FigQ, 
the pod instead is able to overcome this global flow coher- 
ence and to selectively associate the spectral components 
to the adequate spatial coherent structures. This result 
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Data set snaphots number 

FIG. 5: Convergence test of mode 1 passing from p samples 
to p + 1 in the data set. See text for description of the r) 
criterion definition. 



henceforth naturally confirms the mixing layer origin of 
the most energetic spectral component. 

At this step, it might be interesting to briefly dis- 
cuss some critical points of the technique. First, be- 
cause the method is aimed to track out coherent pat- 
terns encountered within a flow (coherent with respect 
to the point- wise correlation matrices), it is important 
for the statistical flow properties to be stationary. As 
a consequence, the data set must possess a sufficiently 
important number of independent realizations so as to 
ensure the convergence of the decomposition towards the 
real pod modes. We have checked that for data set of 
less than 400 samples, the third mode fairly mixes both 
shear layer and cavity structures, resulting in its time 
amplitude fourier spectrum to the occurrence of the 13.5 
Hz-peak — strongly weakened here in mode 3 when us- 
ing 600 samples. Secondly, from an experimental point 
of view, each sample composing the data set should share 
identical (statistical) properties; as a consequence, when 
directly working on instantaneous snapshots of the flow, 
particule feeding should remain homogeneous in time, the 
average intensity and coherent structure resolution being 
modified as the feeding is varying — therefore biasing the 
statistical representativity of the samples pj . There are 
no systematic test to decide whether statistical conver- 
gence has been reached or not. We however plotted in 
FigO the average difference r\ between two modes with 
respect to the number of snapshots contained in the data 
set: 



v(p) = 



W +1 (r)|-|^(r)| 



dr, 



where ip? (r) is the n th -topos computed using p snapshots 
in the data set for the single ^-component of the velocity. 
Note that we had to deal with the absolute value of the 
topos so as to get rid of the sign, since it was observed 
cyclic global sign inversions from vb\ to without 
deep modification of the velocity structure. We see from 
FigElthat convergence is ensured for mode 1 with p ~ 200 
flow realizations. 

The study reported here in fact brings another very in- 
teresting insight from an experimental point of view. It 
indeed shows that, although the velocity field is spatially 



fully 3D and characterized by 3 components [lfj, a 2D 
pod calculation (performed in a plane), over 2 velocity 
components, is able to separate the two intuitive regions 
of interest in the flow (namely the mixing layer and the 
cavity vortices), which therefore strongly simplify any ex- 
perimental protocol, in that a classical PIV (in a plane, 
over 2 velocity components) is sufficient to track out the 
coherent structures and their dynamical features, with- 
out having to call upon 3D PIV techniques, much heavier. 
We have checked that the results were very similar when 
using 1 or 3 velocity components instead of 2. Moreover, 
the 3D calculation of the pod modes confirms all the re- 
sults provided by the 2D analyses; 2D cuts out of the 3D 
modes look very similar to our (intrinsically) 2D modes, 
and their amplitude spectral distribution are comparable 
as well (see [l(j). 

In conclusion, a pod technique has been applied with 
success to discriminate the relevant dynamical features of 
the coherent structures present in the flow over an open 
cavity. The processing time revealed to be of the order 
of 30 s for about 600 samples of size 37300 pixels, and 
grew up to 11 min when applied to about 300 experimen- 
tal PIV samples of size 241.800 pixels (N = 260 x 465). 
However, in most experimental applications, the whole 
field resolution, or the whole picture area, are not re- 
quired to get the expected results, and it is expected 
that the technique could efficiently be applied to a panel 
of other open flows presenting self-sustained oscillations. 

Matlab programs can be obtained from the authors 
upon request. 
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